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ABSTRACT 

The primordial internal structures of gas giant planets are unknown. Often 
giant planets are modeled under the assumption that they are adiabatic, con¬ 
vective, and homogeneously mixed, but this is not necessarily correct. In this 
work, we present the hrst self-consistent calculation of convective transport of 
both heat and material as the planets evolve. We examine how planetary evolu¬ 
tion depends on the initial composition and its distribution, whether the internal 
structure changes with time, and if so, how it affects the evolution. We consider 
various primordial distributions, different compositions, and different mixing ef- 
hciencies and follow the distribution of heavy elements in a Jupiter-mass planet 
as it evolves. We show that a heavy-element core cannot be eroded by con¬ 
vection if there is a sharp compositional change at the core-envelope boundary. 

If the heavy elements are initially distributed within the planet according to 
some compositional gradient, mixing occurs in the outer regions resulting in a 
compositionally homogeneous outer envelope. Mixing of heavy materials that 
are injected in a convective gaseous envelope are found to mix efficiently. Our 
work demonstrates that the primordial internal structure of a giant planet plays 
a substantial role in determining its long-term evolution and that giant plan¬ 
ets can have non-adiabatic interiors. These results emphasize the importance 
of coupling formation, evolution, and internal structure models of giant planets 
self-consistently. 


Introduction 


The existence of a core in giant planets and its physical properties has important con¬ 
sequences for our understanding of planet formation and internal structure (|Helled et ah 


2014). Most investigations of giant planet structure and evolution have assumed that the 


planet consists of several distinct regions which remain chemically homogeneous during the 


long-term evolution (Baraffe et al. 2008; Fortney et al. 2007 Vazan et ah 2013). This is. 


however, a simplifying assumption and in reality it is possible that different regions of the 
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interior will mix or separate as the planet evolves. Changes in composition can also arise 
when one material becomes soluble in another. Recent studies of the solubility of analogous 
phases in the planetary interior have shown that H 2 O ([Wilson fc Militzer||2012a), Si02 (Wil¬ 


son & Militzer 2012b), and iron (Wahl et ah 2013) dissolve in liquid metallic hydrogen at 


high temperature and pressure. These calculations imply that giant planet interiors are in 
thermodynamic disequilibrium with surrounding layers, promoting redistribution of heavy 
elements. 

Alternatively, as the planet evolves parts of the interior may enter the pressure-temperature 
regime where one material is immiscible in another. An example of this would be helium 


rainout during the evolution of a gas giant (Stevenson & Salpeter 1977 Fortney & Hubbard 


2003). This leads to a helium-rich region just outside the core. On the other hand, Steven¬ 


son 


(1982) has suggested that (part of) the planet’s core can be soluble in the surrounding 


gaseous envelope, and that during the formation processes, some of this core material might 
be dredged up and mixed with the surrounding hydrogen and helium. This effect, however. 


is expected to be small (Lissauer & Stevenson 2007). 


The distribution of the heavy elements in the planet, whether concentrated in a central 
core or mixed into a hydrogen-helium envelope, can be of great importance for determining 
planetary thermal evolution and has consequences for the total amount of heavy elements in 


the planet deduced from static models (Stevenson 1985 Chabrier & Baraffe 2007 Leconte 


& Chabrier 2012, 2013). Rearrangement of the heavy elements has important consequences 


for the understanding of the planet formation process. For example, core erosion and mixing 
could explain how massive planets formed by core accretion can have small cores. Such a 
process might also explain, for example, why Jupiter appears to have a smaller core than 


Saturn (Stevenson 1982; Guillot et ah 2004). Furthermore, composition redistribution also 


changes the moment of inertia of the planet, making it harder to deduce the formation 


mechanism even with the more precise measurements expected in the future (Helled 2012). 


An additional assumption that is often made is that there is efficient large-scale convec¬ 
tion in giant planets. However, there are several situations where this assumption is incorrect 

The existence of a core 


Guillot et al. 

1994 

Saumon & Guillot 

2004; 

Fortney et al.| 

2011 


or other gradient in the distribution of heavy elements, in particular, weakens the adiabatic 
assumption. When a stabilizing composition gradient exists, there is no large-scale convec¬ 
tion, and the efficiency of the heat transport can strongly be reduced. Thus, while thermal 
evolution models have been calculated for giant planets, these models generally consider the 
effect of convection on the evolution of the temperature gradient, but not on the evolution 


of the heavy element distribution (Fortney et ah 2011 Chabrier & Baraffe 2007 Mirouh 


et al. 2012). Leconte & Chabrier (2012| 2013) emphasize the importance of a composition 
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gradient on the thermal evolution of the planet, but do not actually follow the evolution 
of the composition gradient. Since convection transfers both heat and material, and since 
the convection itself depends on both the thermal and compositional gradient, these two 
transport processes must be considered simultaneously. 

In this work we present the first self-consistent calculation of convective transport of 
both heat and material. The criterion for the onset of convection is taken to be the Ledoux 
criterion. The Schwarzschild criterion is also considered for comparison. In each case, within 
the convective regions the transport of both heat and material is computed using the mixing 
length recipe. We then calculate full evolutionary models and demonstrate how the compo¬ 
sitional structure changes with time, and how this affects the radius and luminosity of giant 
planets as they evolve. 


2. The method of calculation 

The structure and evolution calculations are carried out using a planetary evolution 
code that solves the equations 
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where V = dlnT/dlnp is the ’temperature gradient’, which may be radiative (and conduc¬ 
tive) or convective, depending on the planet’s local properties (Appendix A). In the last 
equation, Xj is the mass fraction of the j’th species, and 


F, = -{ATrr^pfD^ (6) 

is the corresponding particle flux through the sphere containing mass m. The convective 
diffusion coefficient D (Appendix A) is non-zero in convective regions, otherwise it vanishes. 
The boundary conditions at the center are 


m = r = L = Fj = 0. 


( 7 ) 
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At the surface, which we choose as the photosphere, 

m = M, kpg = {GM/R^)ts, L = Fj = 0, (8) 

where r^, the optical depth of the photosphere, is set to be 1. 

Eqs. 0 -(§ are replaced by difference equations on a grid of n ~ 100 points. In order 
to get an optimal grid resolution during evolution, we add the second-difference equation 


fr+i - fi = fi - fi-u i = 2,...,n-l 


( 9 ) 


where 

/ = {rn/MY^^ — Cl Inp -|- C2 ln(r^ -|- C3), (10) 

and the c’s are suitably chosen, positive, constants. Eq. § forces equal steps of near 
the center, and equal steps of Inp or Inr—whichever step is larger—in the outer region. The 
difference equations are all solved simultaneously for each time step. Further details of the 


code can be found in Kovetz et al. (2009). 


In this work we use the equations of state tables from the work of Saumon et al. (1995) 
for hydrogen and helium. For the high-Z material, which is assumed to be either “ice” 
(H 2 O) or “rock” (Si 02 ), we computed an equation of state based on the quotidian equation 


of state (QEOS) of More et al. (1988), as described in 

Vazan et al. 

(2013). The opacity is 

computed using the radiative opacity tables of Pollack et al. (1985 

for solar-composition. 

and the thermal conductivity tables oflPotekhin et al. 

(1999). 



Unless otherwise stated, the onset of convection in the planetary interior during evolu¬ 
tion is determined by the Ledoux’s criterion. 


AV = — Va — Vx > 0, (11) 

where V_r and Va are the radiative and adiabatic temperature gradients, respectively; and 




d\nT{p,p,X) dXj 
dhip 


dXj 


( 12 ) 


If Vx = 0, convection sets in when V_r > Va, the usual Schwarzschild criterion (?). Fur¬ 
ther details, in particular the calculation of V in eq. (|^, as well as that of D in eq. (§, 


are given in Appendix A. We note that, in accordance with eq. (11), steep enough com¬ 
position gradients—in a direction that renders Vx positive—will cause convection to cease 
completely. If, on the contrary, AV < 0, the actual gradient is the radiative gradient 
and the diffusion coefficient is zero: there is neither convective heat transport, nor any mix¬ 
ing. We do not consider semi-convection—that is, the situation Va < Vr < Va + Vx—or 
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double-layered convection in this study. Thus, our solution (when using Ledoux criterion) 
can provide a lower bound for the mixing efficiency inside the planet. 

Although the mixing is mainly controlled by the initial composition gradients and the 
convection criterion, the rate of mixing can vary with the model parameters. One of the 
most important is the mixing length parameter, a = ^/Hp, which determines the rate at 
which the material is mixed, since it affects the convective diffusion coefficient D. In this 
work we consider different values for a, allowing it to change by several orders of magnitude. 
The planetary evolution and the mixing of heavy elements for different primordial internal 
structure and mixing length parameters are described below. 


3. Heavy-element redistribution 

As we have noted above, the presence of convection and its strength depend on the 
temperature and composition distributions. In order to understand the consequences for 
different primordial internal structures, we consider three representative cases: (1) A planet 
with a core consisting of pure high-Z material, surrounded by an envelope of hydrogen and 
helium in solar ratio. (2) A planet where Z, the mass fraction of heavy elements, decreases 
continuously from a value of Z=1 at the center to Z=0 at the base of the envelope. (3) A 
solar-mix hydrogen-helium planet where the high-Z material is distributed as a Gaussian 
around some point in the interior. This might be taken to represent the injection of material 
by a giant impact in the early history of the planet. Below we present the results for these 
cases. In all the calculations we consider a IMj planet. Unless otherwise noted the high-Z 
material is represented by water and the criterion for convection is Ledoux’s. 


3.1. Core - envelope 

Fig.|I| shows results for three cases of a core-envelope structure. The black dash-dot 
curve represents the initial Z as a function of normalized mass inside the planet, while the 
blue curve shows the Z distribution after 10^° years of evolution. Three different cases are 
presented: The hrst is a massive heavy element core of 0.5 Mj, surrounded by a hydrogen- 
helium envelope (left). The second is a smaller heavy element core of 0.2 Mj , surrounded by 
a hydrogen-helium envelope (middle). The third is a massive heavy element core of 0.5 Mj , 
surrounded by a hydrogen-helium envelope, but here we have mixed an additional 0.01 Mj of 
high-Z material. This material is distributed along a linear gradient starting at Z=1 near 
the core and going to Z=0 at m/M = 0.6. This is meant to represent an initial state that 
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results from some previous heavy-element settling process and/or diffusion. The bottom 
panels show the density (blue) and temperature (green) as a function of normalized mass for 
the different cases after 10^° years. The mixing length parameter in these simulations was 
a = 0.5, but the results are not sensitive to the exact value chosen since a only determines 
the rate at which the material mixes and not the actual onset of mixing. 

In the case of a core-envelope conhguration (left and middle panels), the sharp compo¬ 
sition gradient stabilizes the planet against convection and a radiative layer is formed at the 
core-envelope boundary. As a result, although large parts of the envelope do convect, con¬ 
vection is inhibited near the core, and the distribution of heavy elements within the planet 
remains essentially unchanged. Even after 10^° years of evolution the outer envelope heavy 
element mass fraction is negligible. A consequence of this inhibition of convection is that 
the planet retains its heat longer. Therefore, the memory of the initial temperature prohle, 
which is generally lost after ~ 10'^ years, can persist for as long as ~ 10^° years. 

Lowering the mass of the core does not cause it to mix more readily. This is demon¬ 
strated in the middle panel, where even a 0.2 Mj core remains unmixed after 10^° years. 
However, second order processes, like diffusion or chemical solubility, might still slowly erode 
the core and establish a composition gradient in the envelope. In the right-hand panel of 
Fig.|Il we simulate a case where a thin region with a composition gradient exists around 
the core. In this case, as shown in the hgure, the outer part of this region becomes convec¬ 
tive during evolution, and the high-Z material mixes uniformly through the outer envelope. 
This mixing, however, occurs only in the outermost part of the planet , where the thermal 
expansion coefficient is high, and does not affect the core directly. It therefore appears that 
core erosion cannot proceed via rapid upward mixing when there is a sharp compositional 
gradient, but only through much slower diffusive processes. 


3.2. Gradual composition gradient 

Once the composition gradient is reduced from being a sharp step to something more 
gradual, convective mixing may come into play in different regions in the planet. A typical 
case is shown in Eg. The high-Z distribution is shown for different times in the evolutionary 
sequence. The color of the curves shows the associated temperature profile. 

As can be seen from the hgure, the upper 20% of the mass is nearly completely mixed 
after 10® years, while the innermost 80% or so maintains its original Z-prohle. There is 
an intermediate narrow region between the radiative and convective where convection is 
not completely suppressed. The steps in heavy element composition in this region may be a 
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numerical artifact. In any case their effect on the subsequent structure and evolution is small, 
as we show below (section 3.2.2). In the outer regions convection becomes strong enough 
to completely mix the overlying material, and indeed the heat transport tracks the material 
transport closely. The inner 80% of the planet stays at roughly the same temperature while 
the outer 20% cools as the mixing proceeds. Fig. shows the effect of changing the slope 
of the composition gradient. For all three panels Z=1 at the center and decreases linearly 
with radius to a surface value determined by the desired total Z for the planet. In all cases 
the original heavy element distribution is retained in the inner region of the planet while the 
outer part becomes homogeneously mixed. 


3.2.1. Dependence on mixing length parameter 


The mixing length parameter used in Fig. and is a = 0.5. However, since the 
mixing rate depends on the value of a, it is important to investigate how the mixing of the 
heavy elements in the planet changes for different values of a. In stellar evolution models 
the mixing length is typically taken to be of the order of the pressure scale height, so that 
a ~ 1. In planetary interior conditions, however, the mixing parameter could be signihcantly 
smaller (see, e.g., Leconte & Chabrier 2012). Fig. shows how the choice of a affects the 


convective mixing during the planetary evolution. For a = 0.5 (upper left) mixing is efficient 
and remains so for a > 10“^. Lower values reduce the mixing efficiency until for a = 10“® the 
internal compositional structure is essentially unchanged even after 10^*^ years of evolution. 


For such low a values (Leconte & Chabrier 2012) radiative transfer is more effective than 
convection. 


3.2.2. Dependence on number of mesh points 

In order to assess the numerical stability of the results we ran the model of Fig. |^with 
150 and 500 mesh points for 10^° years. The results are shown in Fig. As the resolution 
is increased (i.e., more mesh points) the number of steps in the high-Z material increases in 
the region m/M = 0.6 — 0.8, although the position of these steps remains roughly the same. 
The location and extent of the radiative and convective regions are basically unchanged. 
As a result the temperature prohle for both cases is essentially the same (lower panel of 
Fig. 1^). The large step at the outer edge of the planet is slightly reduced. This is because 
the increased resolution allows us to £x the outermost convective region more precisely. The 
actual difference is relatively small, however, and the value of Z in this outer convective 
region decreases from 0.015 to 0.010. As a result, the radius as a function of time remains 








essentially unchanged as well. We can therefore conclude that the sensitivity of the calculated 
evolution to the number of mesh points is relatively small. 


3.2.3. Dependence on composition 


Heavy elements in giant planet interiors can be mainly icy, or rocky. In reality, planetary 
interiors probably contain both. Necessarily, different thermodynamic properties of each 
material will affect the convective mixing. In addition to the known effect that the radius 


varies depending on whether rock or ice is assumed for the heavy material (Vazan et ah 


2013), there can be an effect on the evolution via mixing. For this reason, we compare the 


evolution of a IMj planet with an icy to one with a rocky core. Fig. shows the difference 
in mixing when using H 2 O vs. Si02 for a similar initial model. Here too we use a = 0.5 
and start with a more gradual gradient in heavy elements. The initial high-Z distribution is 
shown by the dotted-black curve, and the distribution after 10^° years is shown by the blue 
and red curves for H 2 O and Si02, respectively. As expected, mixing for the case of Si02 
is less efficient than for H 2 O because Si02 has a higher molecular weight. For the case of 
H 2 O the outer 50% of the planet (by mass) is mixed and Z=0.08 in this region. For the 
Si02 case only the outer 30% of the planet is mixed and Z=0.03. In principle, the increased 
metallicity of the envelope due to mixing can also affect the planetary opacity, which could 


in turn affect the subsequent evolution of the planet and its observable parameters (Burrows 


et al. 2007; Vazan et ah 2013). We hope to address this topic in future work. 


3.3. Gaussian distribution 

The third set of cases we explore is a Gaussian distribution of high-Z material centered 
around some mass level in the planetary hydrogen-helium envelope during the early evolution. 
The results for this configuration for different masses and locations of high-Z material is 
presented in Fig. The left-hand panel presents the evolution for lOMQof H 2 O injected 
at m/M = 0.3. The dash-dot curve is the initial distribution. Mixing is found to be both 
downward and upward, and after 10^ years there is already considerable mixing (blue curve). 
After 10® years, the envelope is thoroughly mixed (red curve). 

Injecting 16 M 0 (middle panel) and 20 M 0 (left panel) gives qualitatively similar results. 
In addition, we find that the location in which the material is injected does not affect the 
final results of complete mixing. The effect of the composition gradient on convection for 
these cases is small. We compared the evolution using the Ledoux criterion with that using 
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the Schwarzschild criterion (not shown) and found very similar behavior. The low fraction 
of high-Z material (lower than Z=0.3) expected from the dissolution of even a 16 M 0 rocky 
or icy object has hardly any effect on the convection criterion. The mixing has no real 
(observable) effect on the physical and thermal parameters of the planet during evolution, 
but it does affect the planetary envelope enrichment. The evolution of radius, temperature, 
density prohles, as well as the planetary luminosity, remains the same whether the high-Z 
material is originally completely mixed, or is in a Gaussian distribution around some mass 
point. 


3.4. Radius and luminosity evolution 

Since both the radius and luminosity of exoplanets can be measured, and they play key 
role in the characterization of the planets, we next investigate how the internal structure 
and mixing affect these properties. The change of the planetary radius (top) and luminosity 
(bottom) with time for a continuous compositional gradient is presented in Fig. The radius 
is given in units of Jupiter radius (Rj) and the luminosity in units of solar luminosity, (Ls^^). 
We compare two cases: the blue curves show the evolution when the criterion for convection 
is taken to be the Schwarzschild criterion while the red curves show the evolution using 
the Ledoux criterion. For each of these two cases we further consider only heat transport, 
that is, convection with zero diffusion coefficient D (solid curves), and heat together with 
high-Z material transport (dotted curves). As expected, the Ledoux criterion leads to slower 
contraction due to less efficient heat transport, and results in a larger radius. When material 
mixing is allowed {D 7 ^ 0) the planetary radius decreases during the evolution, mainly due 
to the increased density of the envelope caused by the enrichment in heavy elements. The 
difference in luminosity between the various cases is about a factor of 2-3 throughout the 
evolution. As expected, giant planets that are modeled using the Schwarzschild criterion 
are found to be more luminous than planets modeled on the Ledoux criterion, due to more 
efficient heat release. Since mixing occurs only at a later stage the luminosity (and radius) 
for the cases with and without mixing for both Schwarzschild and Ledoux is the same up 
to about 10 ^'® years, from that point on the difference in luminosity (and radius) becomes 
more apparent. For both the Ledoux and the Schwarzschild scenarios the mixing reduces the 
radius, leading to a lower luminosity. However, while the luminosity for the Schwarzschild 
scenario decreases when mixing is included, for the Ledoux case it increases. This means 
that the effect of the slower heat release in the case of Ledoux is more important than the 
reduction in radius (surface area) due to the heavier envelope. 

The structures of core-envelope and the completely mixed planet are often considered 
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as the extreme cases which bracket intermediate structures. In Fig. [^the radius (top) and 
luminosity (bottom) evolution of planets of IMj , with the same initial Z fraction (Z=0.35), 
but different initial distributions, is presented. We consider three different cases: (1) a 
planet with a continuous composition gradient (red-dash), (2) an initially fully mixed planet 
(green), and (3) a pure high-Z core -|- hydrogen-helium envelope (blue). The fully mixed 
planet has large convective regions during its evolution (no compositional gradients), and 
heat can be rapidly transported out of the planet. For the core-envelope case convection 
is suppressed at the core envelope boundary, so that the core cools more slowly, but the 
envelope convects and loses heat fairly rapidly. Its radius stays larger than the fully mixed 
case, in part for that reason. However, the planet with the gradual gradient in composition 
has convection suppressed throughout most of the volume, and therefore cools more slowly 
than either of the previous cases. As a result, its radius is larger that either of the two 
’bracketing’ cases, bottom panel of the hgure. The difference in planetary luminosity for the 
different cases, as shown in the bottom panel of the hgure, are related to the heat release rate 
of the planets. The homogeneous planet, which has the most efficient heat release, has the 
higher luminosity during most of the evolution time. The gradual planet, which releases heat 
less efficiently, has the lower luminosity. The correlation changes after a few billion years, 
since the difference in radius also affects luminosity: the fully-convective homogeneous case 
has a smaller radius. 


The influence of mixing on the energy budget of the planet is negligible. This is because 
the amount of material that mixes is relatively small, so that the amount of energy taken 
up by raising the high-Z material out of the potential well does not affect the energy of the 
planet. One might expect that the upward mixing of high-Z material would convert thermal 
energy into gravitational energy thus temporarily speeding up the cooling of the planet. 
Such an effect, however, if present, is small, and there is hardly any noticeable difference in 
the thermal energy between the mixed and unmixed case. In fact, the mixing causes the 
total radius of the planet to decrease more quickly, as can be seen in Fig. and so the 
gravitational energy becomes more negative for the case with mixing. 


Finally, Fig. 10 shows the evolution of IMj planet with different mixing length param¬ 
eters. Lower a values result in slower convection, which leads to slower cooling and larger 
radii but lower luminosities. Since the mixing length parameter, a, affects the convection 
and mixing velocity, but not the occurrence of convection, different internal structures are 
developed during the evolution, leading to different convective regions in each model. This 
results in different cooling histories and radius and luminosity evolution. Thus, there is 
no simple correlation between a and the planetary radius/luminosity, as can be seen from 

Fig.IT^ 



11 


It is clear from the figures that the radius and luminosity at a given time for a planet do 
not depend only on its composition, but are also dependent on the distribution of the heavy 
elements as well as on whether mixing of material occurs as the planet evolves. Therefore, 
one must be careful when using the simple mass-radius diagram for characterizing giant 
planets and inferring their bulk compositions. 


4. Discussion and conclusions 


We investigate how convective transport of heat and high-Z material affects the evolution 
of giant planets. The evolution is followed allowing for both heat transfer and mixing of high- 
Z material. We consider various initial distributions of heavy elements, and hnd that the 
primordial internal structure has an important role in determining the subsequent evolution 
and hnal internal structure. We use both the Schwarzschild and the Ledoux criteria for 
the onset of convection. When the Ledoux criterion is used, we hnd that, in the case of a 
distinct core-envelope boundary, core erosion is not possible due to the sharp transition in 
composition, which suppresses convection. 

An initial interior structure in which the heavy elements are distributed with a gradual 
composition gradient within the planet can become homogeneously mixed in the outer regions 
due to convection. The amount of mixed material depends on the mixing length parameter 
as well as on the composition of the heavy elements. The main effect of convective mixing 
is to enrich the envelope in high-Z material. This reduces the radius of the planet during 
the long-term evolution. It is often assumed that the evolution of a planet can be bracketed 
by two end-member models: a planet where the high-Z material is fully mixed and a planet 
where the high-Z material is completely in a core. We have shown that a planet with a Z 
distribution that varies continuously within the interior can have a radius larger than both 
the end-members due to less efficient heat release. Because the details of the mixing depend 
on the initial composition and its depth dependence, and because these details influence 
the evolution and the hnal structure of the planet, we conclude that the primordial internal 
structure of giant planets, in particular the initial high-Z mass distribution, is important and 
can affect the subsequent evolution. 


The primordial internal structure of a giant planet depends on its formation mechanism. 
Although both core accretion and disk instability can lead to various compositions and 
internal structures ( [Helled et ah 2014), it is often felt that giant planets formed by core 
accretion will have relatively massive heavy-element cores, while giant planets in the disk 
instability scenario are expected to be more homogeneously mixed. In a way, the core- 
envelope internal structure we consider could be a representative of core accretion while the 







12 


gradual change in Z could represent disk instability planets if efficient planetesimal accretion 
takes place. However, a gradual distribution of heavy elements is possible in planets formed 


by core accretion, depending on where the heavy elements are deposited (laroslavitz & 


Podolak[2007) and whether major mixing occurs during runaway gas accretion. It is therefore 


clear that we need to establish a better understanding of the expected composition and heavy 
element distribution from giant planet formation models. Since the primordial distribution 
of heavy elements depends on the formation process, and since the evolution depends on the 
primordial internal structure, and at the same time, the internal structure can change with 
time due to mixing and settling, a coherent picture of giant planets can be achieved only by 
linking the three aspects of origin, evolution, and internal structure self-consistently. 
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Appendix A: Ledoux’s criterion and the mixing-length recipe 


Ledoux’s criterion for the convective stability of a hydrostatic configuration is usually 
stated for mixtures in which the equation of state (EOS) depends on composition through 
the molecular weight /r. In an ideal gas, for example, the EOS is p = {RG/fi)pT , where Rg 
is the universal gas constant. In other cases, for example when the EOS for a mixture with 
mass fractions (Ai, is 

1 = 

P ^ PjiP,T) 

(the additive-volume law), Ledoux’s criterion requires some generalization. We prefer to 
derive it ab initio. 


The criterion for convective stability is 

dR ^ dp{p,s,X) 
dp ~ dp 


m 


where the derivative on the right hand side is taken at constant specific entropy s and 
constant mass ratios Xj (collectively denoted by X). On the left-hand side, we have the 
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slope of the configuration’s (p, p) structure line. Now 


ip = 


dp 


dT 


dX. 


{A2) 


where the repeated index in the last term implies summation over all specie^ Thus, the 
left-hand side of (Al) is 


dp _ dp{p, T, X) dp{p, T, X) dT dp{p, T, X) dXj 
dp dp dT dp dXj dp 


According to (A2), 


dT{p,p,X) dp{p,T,X)/dX, 


dXi 


dp{p,T,X)/dT- 

If we replace p{p,s,X) by p{p,T,X), where T = T{p,s,X), the right-hand side of (Al) 

bGCOIHGS 

dp{p,s,X) dp{p,T,X) , dp{p,T,X)dT{p,s,X) 


dp 

With the standard notations 


dp 


V = 


dlnT 


+ 




dlnp ’ 

as well as the additional notation 

d\nT{p,p,X) dXj 


dT dp 

d\nT{p,s,X) 
dlnp ’ 




dX, 


dlnp’ 


(A3) 


the criterion (Al) for convective stability becomes 

dp{p,T,X) 


dlnT 


(V — Va — Vx) > 0. 


(A4) 


For materials satisfying Q = —dinp{p,T,X)/dlnT > 0 (positive volume expansion - rather 
than contraction - coefficient) this becomes 


V — Va — Vx < 0, 


(A5) 


which is the desired, general form of Ledoux’s criterion for convective stability. In a layer of 
uniform composition Vx = 0, and (A5) reduces to Schwarzschild’s criterion V — Va < 0 . 


^In differentiating with respect to the Xj’s, we disregard the fact that = 1- This is guaranteed by 

eqs. (5)-(6): since the diffusion coefficient D is common to all species, we have ^ whenever ^ 

throughout the planet. According to eq. (5), then, dY^^j/dt = 0 everywhere, and Y^i — ^ ^ih continue 
to be satisfied throughout. The requirement that Y W “ ^ everywhere is in fact an initial condition. 
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If, as in an ideal gas, the EOS depends on the mass fractions only through the molecular 
weight, then 

c?hiT(p, p, X) dXj 91nT(p,p, p) 5/i dXj 
dXj d\np dp dXjdlnp^ 

and the inequality (AS) assumes its familiar form 


V - ^ < 0. (A6) 

a Inp 

since, in a perfect gas, d\nT{p,p, p)/d\np = 1. 

Linear stability analysis predicts that, when the inequality (Al) is violated, small per¬ 
turbations will grow exponentially with time (?), but the actual amplitudes reached lie 
outside the linear analysis: they are limited by the non-linear terms that were left out by 
that analysis. 

If the total luminosity is due to radiation (and conduction). 


L 




c dpR 
Kp dr ’ 



The ‘temperature gradient’ dlnT/dlnp that results from this equation and the hydrostatic 
equation ([^ is 

dlnT p kL , 

(A7) 


Vk = 


p 


kL 


d\np Apn AtvcGui 
Now, if this temperature gradient violates the stability criterion, that is, if 


Vp — Va — Vx > 0, 


then convection will set in: a part of the luminosity will then be carried by a convective flux, 
and the rest by radiation/conduction, requiring a V which is less than Vr. 


For any V > Va + Vx the excess in dp/dlnp of a convective element, relative to its 
suroundings, is given by the left-hand side of (A4). The density excess after travelling a 
‘mixing-length’ distance i is obtained by multiplying by id\np/dr. The resulting density 
excess is therefore 


5p ~ 


dp{p,T,X) i 
dlnT Hr, 


where Hp is the pressure scale height, dehned by ^ 
temperature excess then follows from eq. (A 2 ): 


(V- Va-V x), 

= —dlnp/dr. 


The corresponding 


6\nT ~ 


lip 
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Actually, however, convective elements do not move adiabatically, since they continually lose 
energy to their surroundings by radiation or conduction. We therefore replace in the 
foregoing formulae by the (as yet unknown) element’s temperature gradient V'; 




(A8) 


hlnT~ -V). 
tip 

From (v ■ grad)v ~ g(51np we obtain an estimate for the convective speed. 


(yl9) 


~ gU\np ~ gH^Q( — )\v -V- V x] 

11 r) 


(9110) 


The convective energy flux is Fc = pv6h, where 6h is the excess heat per unit mass, given 
by the change in the specific enthalpy h = u + p/p, where u is the specihc internal energy. 
At a given pressure, we have 


dh = 




dT 


dX, 


^ , dh{T,p,X)^^ 
= CpT d In T H-—- dXj , 

(y d 


where, again, the repeated index implies summation over all species. Thus, according to 
eqs. (A3) and (A9), the excess specihc heat is 


6h ~ CpT—(V-V') + 

til 


i dh{T,p,X) dXj 


dX, 


^ c,rT(v-v'-Vx) + 9y(c,r 

lip lip 


d\np 

d\nT{p,p,X) dh{T,p,X)\ dXj 
dXj dXj ) dinp 


In the last sum over compositions we expect some cancellation: after all, since '^Xj = 1, 
the dXj/dlnp^s cannot all have the same sign. Moreover, in a mixture of perfect gases the 
sum is identically zero because each one of the coefficients of the dXj/dlnp^s vanishes! In 
what follows we shall neglect this sum and estimate the excess heat per unit mass by 


^CpT(V-V'-Vx). (All) 

lip 

The resulting convective energy hux Fc = pvSh is then (cf. eq. (AlO)) 

Fc ^ {gHpQY/^pCpT{^)\x - V' - (A12) 

lip 
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Two equations are needed in order to determine the unknowns V' and V. We first con¬ 
sider the ratio of excess energy of an element after travelling a distance ^ relative to its 
surroundings, which—according to eq. (All)— is proportional to (V — V' — Vx), to the 
energy lost by radiation/conduction during the time of travel ^/v^ which is proportional to 
(V — Va — Vx) — (V — V' — Vx)- An estimate (Mihalas 1978) of this ratio yields the first 
equation, 

V' - Va 16ctT4 


(V - V' - Vx)iA pCpT{gH^QYI^{l/H,yrH' 
where a is the Stefan-Boltzmann constant and th = KpHp. The second equation is 


(^13) 


(AU) 

where Fr = —c/ {Kp)dpji/dr = 16(jT^V/(Str) is the radiative flux, and F = 16(jT^Vr/(Str) 
is the total flux. A mixing-length recipe consists in replacing each ~ in the foregoing relations 
by an equality, followed by a numerical factor. According to Mihalas (for example) 

V = {gHpQ/8Y/\i/Hp){V -V - 
F, = {gHpQ/2>2f/^pCpT{l/Hpfiy-W'-Wxf\ 

V' - Va _ 32y2aT" 

(V - V' - Vx)V2 pCpT{gHpQYIYllHpYTH ~ ' 

Eq. (A14) and the last two equations lead to the cubic equation 

2h'x + = 1, (A15) 

46' 

where 

b' = 6/(Vr - Va - Vx)'/^ = (V - V' - Vx)/(Vr - Va - Vx). 

The cubic equation is the same as the one used in the Schwarzschild context, but b' and 
X are differently defined (except when Vx = 0). Finally, the solution of (A15) yields the 
required temperature gradient to be used in eq. (|^; 

V = Va + Vx + (Vr — Va — Vx)(26'a; -1- xY- (A16) 


The last equation, again, reduces to the one used in the context of Schwarzschild’s criterion 
when Vx = 0. Efficient convection (6' << 1) leads to V ~ Va + Vx; when convection is 
inefficient, V ~ Vr. The ratio of convective to total energy flux is 


Vr — V Vr —Va —Vx 3 


-x^. (A17) 

Vr Vr 46' ^ ^ 

Finally, mixing in a convective layer is effected by eq. (|^, where the ‘convective diffusion 

coefficient’ D —common to all species—is a fraction of vi (as if the convective elements 

were gas particles, moving at a mean speed v with mean free path i). In this study we set 

D = 0.1u£ = 0.1{i/Hp)vHp. In a radiative layer V = Vr and D = 0. 
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Fig. 1.— Heavy element (H 2 O) mass fraction vs. normalized planetary mass for a 1 Mj planet 
after 10^° years of evolution for a planet consisting of a higli-Z core and a hydrogen-helium 
envelope. Initial internal structure (dashed dotted) and after 10^° years of evolution (solid) 
are presented. Left: Me = 0.5 Mj; Center: Me = 0.2 Mj; Right: Me = 0.5 Mj with a 
small Z-gradient above the core (0.01 Mj). The bottom panels show the density (blue) and 
temperature (green) vs. normalized mass for the different cases after 10^° years. 
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log time= 7 log tlme= 8 



m/M 



logtlme=9 logtlme=10 




Fig. 2.— Evolution of heavy element (H 2 O) distribution for IMj planet in time for a con¬ 
tinuous Z-gradient at different times in the evolution. The step structure (see text) can be 
seen developing. The colors on each curve show the evolution of the temperature prohle. 
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Fig. 3.— Upper panels: Z as a function of m/M after 10^° years of evolution of a 1 Mj planet 
for different initial Z-gradients and total Z with the Ledoux criterion. Shown are cases 
with Z=0.65 (left), Z=0.35 (center), and Z=0.10 (right). Lower panels: Density (blue) and 
temperature (green) after ~ 10^° years for the cases in the upper panels. 
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Fig. 4.— Z redistribution as a function of normalized mass after 10^*^ years of evolution for 
the different values of mixing length parameter, a. 
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150 layers 



500 layers 
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Fig. 5.— Upper panel: Z as a function of normalized mass for a 1 Mj planet 150 mass mesh 
points (left) and 500 mesh points (right). The initial structure is the same for both cases 
and is given by the dashed-dotted black curve. The red curves show Z after 10^° years for 
the two cases. Lower panel: temperature prohle after 10^° yr of evolution, for the above 
planets with 150 (blue) and 500 mesh points. The initial temperature prohle is presented 
also (dashed-dotted). 
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HO SiO, 

2 2 




Fig. 6.— Z as a function of normalized mass after 10^° years of evolution for a IMj planet 
with a total Z=0.35. The initial distribution is given by the dash-dot curve. The final 
distribution is shown for H 2 O (blue) and for Si 02 (red). 
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Fig. 7.— Z as a function of normalized mass for a 1 Mj planet for different Gaussian dis¬ 
tributions. Left: 10 M 0 centered at m/M = 0.3; Center: 16 M 0 centered at m/M = 0.8 ; 
Right: 20 M® centered around m/M = 0.3. The dashed-dotted black curves correspond to 
the initial distribution of the heavy elements. Also shown are distributions after 10^ years 
(blue), 5 X 10^ years (green) and 10® years (red). 
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Radius evolution - 1 Mj planet, gradual, Z=0.35 




Fig. 8.— Radius (top) and luminosity (bottom) as a function of time for a 1 Mj planet, with 
Z=0.35 and an initial continuous Z-gradient (as in fig. center). Shown are evolution when 
using the Schwarzschild convection criterion (dotted blue) and the Ledoux criterion (dotted 
red). The solid curves show the evolution when there is only energy but no high-Z transport. 
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Radius evolution - 1 Mj planet, Zfrac=0.35 



Luminosity - 1 Mj planet, Zfrac=0.35 



Fig. 9.— Radius (top) and luminosity (bottom) as a function of time for a 1 Mj planet for a 
continuous Z-gradient with Z mixing (red dash), pure high-Z core plus solar envelope (blue), 
and homogeneously mixed (green). In all cases Z=0.35. 















Radius evolution - 1 Mj planet, Z=0.35 of HgO, gradual initial 



Luminosity - 1 Mj planet, Z=0.35 of H^O, gradual initial 



Fig. 10.— Radius (top) and luminosity (bottom) as a function of time for a 1 Mj planet with 
Z=0.35 and an initial continuous Z-gradient. Different colors are for different mixing length 
parameters, a. The internal high-Z distribution for those cases is presented in £g. 
















